// ASR trager dan origineel

//----------------------------------------------------------------------------
// M A R I T I M E  R E S E A R C H  I N S T I T U T E  N E T H E R L A N D S
//----------------------------------------------------------------------------
// Copyright (C) 2008 - MARIN - All rights reserved - http://www.marin.nl
//----------------------------------------------------------------------------
// Program    : mermaid
// Module     : vbm
// File       : RrbElwinSolver.h
//----------------------------------------------------------------------------

#ifndef RRBELWINSOLVER_H
#define RRBELWINSOLVER_H

#include "Solver.h"
#include "Array2D.h"

#ifdef EWLAPACK
#include <lapackpp.h> // Lapack++ library // see NopDefSolver.cpp
#endif

#define RRBSOLVER_BASE Solver

/// recursive red black solver by Auke van der Ploeg

class RrbElwinSolver : public RRBSOLVER_BASE
{
    typedef RRBSOLVER_BASE base;

public:

    RrbElwinSolver(Array2D<REAL> *pcc,
              Array2D<REAL> *pcs,
              Array2D<REAL> *pcw,
              int Nx1,
              int Nx2,
              REAL tolerance,
              int maxLevel);

    virtual ~RrbElwinSolver();

    virtual void solve(Array2D<REAL> *pB,
                       Array2D<REAL> *pX);

private:

    /// @cond nevercall
    RrbElwinSolver(const RrbElwinSolver& m);
    RrbElwinSolver & operator =(const RrbElwinSolver& m);
    /// @endcond

    void preconditioner();

    void precsolve(Array2D<REAL> *pz);

    void matvec(const Array2D<REAL> &p,
                Array2D<REAL> *pq,
                Array2D<REAL> *py);

    Array2D<REAL> *m_pcc; /// pointer to MatC
    Array2D<REAL> *m_pcs; /// pointer to MatS
    Array2D<REAL> *m_pcw; /// pointer to MatW
    Array2D<REAL> m_ccc; /// copyMatC
    Array2D<REAL> m_css; /// copyMatS
    Array2D<REAL> m_cww; /// copyMatW
    Array2D<REAL> m_cee; /// copyMatE
    Array2D<REAL> m_cnn; /// copyMatN
    Array2D<REAL> m_Chol; // Cholesky decomposition
#ifdef EWLAPACK
    LaSymmBandMatDouble m_LaChol; // Cholesky decomposition
#endif
#ifdef EWRES
    Array2D<REAL> m_Sc; // original center stencil
    Array2D<REAL> m_Ss; // original south stencil
    Array2D<REAL> m_Sw; // original west stencil
#endif

    int m_dimS; // dimension of last Schur complement
    int m_bandwith; // bandwith of last Schur complement
    int m_nivmax; // maximum level of RRB-iteration
};

#endif // RRBELWINSOLVER_H
